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ABSTRACT 

Inspired by Barnard 68, a Bok globule, that undergoes stable oscillations, we perform 
multi-phase hydrodynamic simulations to analyze the stability of Bok globules. We 
show that a high-density soft molecular core, with an adiabatic index 7 = 0.7 embed¬ 
ded in a warm isothermal diffuse gas, must have a small density gradient to retain the 
stability. Despite being stable, the molecular core can still collapse spontaneously as it 
will relax to develop a sufficiently large density gradient after tens of oscillations, or a 
few 10^ years. However, during its relaxation, the core may abruptly collapse triggered 
by the impingement of small-amplitude, long-wavelength (^6 — 36 pc) sound waves 
in the warm gas. This triggered collapse mechanism is similar to a sonoluminescence 
phenomenon, where underwater ultrasounds can drive air bubble coalescence. The col¬ 
lapse configuration is found to be different from both inside-out and outside-in models 
of low-mass star formation; nonetheless the mass flux is close to the prediction of the 
inside-out model. The condition and the efficiency for this core collapse mechanism 
are identified. Generally speaking, a broad-band resonance condition must be met, 
where the core oscillation frequency and the wave frequency should match each other 
within a factor of several. A consequence of our findings predicts the possibility of 
propagating low-mass star formation, for which collapse of cores, within a mass range 
short of one order of magnitude, takes place sequentially tracing the wave front across 
a region of few tens of pc over 10^ years. 
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1 INTRODUCTION 

Interactions of molecular clouds in the in- 

terstellar medi um _ (I SM) with shock waves 

lStone_&_^Jorma3 Il992l: iK Iein. McKee fc Colellal 1 19941 : 
IXu fc Stond Il9 ^: iNak amura et al.l 2 0061) and tur¬ 
bulence _ (IVazauez-Semadeni. Gazo l fc Scalol 200C; 

Mac Low fc KlessenI l2004l : iLi fc Nakamural 200f ; 

Nakamura &: Lil l2007l) have been extensively studied 
in the past decades. The general consensus from these 
studies is that clouds are most likely to be destroyed 
by non-radiative shocks due to Kelvin-Helmholtz and 
Rayleigh-Taylor instabilities, but can be severely com¬ 
pressed leading to gravitational collapse by weak radiative 
shocks; moreover turbulence generated by proto-stellar 
winds is capable of self-regulation for star formation. While 
these results are valid in highly active regions of ISM, such 
as in giant molecular clouds where nearby proto-stellar 


* E-mail: chiuehth@phys.ntu.edu.tw 


winds, supernova explosion and ionization bubbles are at 
work, the majority part of ISM may on the other hand 
experience only mild disturbances originated from distant 
active regions. Here small starless cores of few solar masses 
in relatively quiescent environments are likely much more 
abundant in the ISM to contribute to the low-mass end 
of the stellar initial mass function. In particular, B68, a 
Bok globule, has been identified to undergo osc illations 
llLada et al.ll2003l : iRedman. Keto fc Rawlin^l2006l ). and it 
is for this population of perturbed low-mass cores that the 
present work aims to focus. 


Theoretical studies of low-mass star formation in qui¬ 
escent environments have historically been analyzed with 
two different model configurations — the ou tside-in collapse 
llLarsonI Il969l : PenstonI 1 19691 : iHunterlIl977l) and the inside- 
out collapse (IShul 1977h . The outside-in collapse model con¬ 
siders an initial uniform cloud that just exceeds the Jeans 
mass and the collapse is initiated near the cloud bound¬ 
ary where the collapse front propagates from outside to¬ 
ward the center. The inside-out collapse model on the con- 

































2 Ui-Han Zhang, Hsi-Yu Schive and Tzihong Chiueh 


trary considers an initial cloud already relaxed to a singular 
isothermal sphere configuration and the collapse is triggered 
at the center where the collapse front subsequently propa¬ 
gates outward much like an avalanche. Weak disturbances 
in ISM were never considered to be important to the pri¬ 
mary physical mechanisms those works put forth. In this 
work, we shall instead demonstrate that a small-amplitude, 
long-wavelength sound may trigger the collapse of an oth¬ 
erwise stable core. This new mechanism involves two key 
components: a soft equation of state in the core and the 
sonoluminescence-type resonant sound absorption. 

It is well known that in a molecular cloud containing 
dust and metal, cooling is dominated primarily by molecular 
line emissions and heating provided by cosmic ray bombard¬ 
ments. While the latter depends linearly on the hydrogen 
density, the former depends more strongly and complicat- 
edly on the hydrogen density. The different density depen¬ 
dence of cooling and heating yields a peculiar equation of 
state, with a drop in the adiabatic index 7 below 1 for a 
hydrogen number de nsity exceeding several times W^/crrd 
dSnaans fc Silkl Bood l . Observation of starless clumps also 
supp orts a decreasing temperature toward the high-density 
core IIBachiller. Guilloteau fc Kahanall987l ). Therefore the 
molecular core must have a rather soft equation of state. 

On the other hand, sonoluminescence originally refers 
to underwater air bubbles that undergo radical compression 
to emit optical lights when subject to the i mpingement of 
sounds (IBrenner, Hilgenfeldt fc Lohsd 12002 1 . When a pla¬ 
nar sound wave in water encounters a soft air bubble, it 
can become a spherical wave propagating into the bubble 
with an increasing amplitude. Mismatch in sound speeds 
may however produce wave reflection, but when the sound 
is in resonance with the air bubble, the wave can be ab¬ 
sorbed with high efficiency. It is in this general context that 
we report this mechanism to be possibly operative in the 
ISM molecular cores. 

We describe our simulation details in Sec. Q- The sta¬ 
bility condition of a soft molecular core is derived in Sec. (Ell- 
Sec. Q presents our main results, addressing the interaction 
of a small-amplitude wave and a single core. Simulation on 
the propagating star formation is demonstrated in Sec. 
where a single propagating wave encounters several cores. 
We give conclusions in Sec. 


2 NUMERICAL SCHEME, EQUATION OF 
STATE AND INITIAL CONDITION 

We demonstrate the core collapse mechanism numeri¬ 
cally. Our simulations are conducted with the GAMER 
code, a graphic-process-unit (GPU) accelerated hydro¬ 
dynamics adaptive-mesh-refinement (AMR) code, which 
outperforms the CPU cod e by a factor of several 
tens in computation speed (Schive. Tsai fc Chiuehl I 2 OIOI : 
ISchive. Zhang fc ChiuelJ 12012 ). The GAMER code solves 
ideal hydrodynamics with self-gravity in conservation forms, 
and it has installed several solvers, such as CTU (corner 
transport upwind), MHM (MUSCL-Hancock) and RTVD 
(Relaxing Total Variation Diminishing). The RTVD solver 
is chosen for this problem. The scheme adopts an approxi¬ 
mate measure for estimating the information speed, making 
it convenient to handle an unusual equation of state adopted 


in this work. With GAMER, we may search for the nonlinear 
stability boundary with hundreds of 3D simulations, much 
like Monte-Carlo simulations, in a sizable simulation box. 
We are also able to achieve (3.145728 x 10®)® effective reso¬ 
lution when the core collapses. 

A simplified initial bare core is adopted, which is to be 
relaxed to an appropriate configuration. This initial dense 
H 2 core has 0.2 pc diameter, ~ 20 K temperature and ~ 10 "^ 
cm“® atomic hydrogen density, surrounded by a voluminous 
HII gas of ~ 10 "^ K temperature and ~ 10 cm“® density and 
separated by an infinitely thin transition layer. 

Our simulation domain is a 37 pc cube with 0.048 pc 
base grid resolution. The base grid resolution is defined as 
level 0 and the n-th refinement is level n, with level n twice 
higher resolution than level n — 1. The space position and 
the density are used as the refinement criteria. For the space 
position, any grid with distance from the core center smaller 
than 2®“"Rc is refined to level n, where n = 1 - 4 and Rc 
is the initial core radius. This condition guarantees the core 
is embedded in grids of sufficiently high resolution. On the 
other hand, the grid with density higher than 37.5 * 4"“^pg 
is also refined to level n, where n = 1 - 12 and pg is HII mass 
density. This condition aims to capture the collapse center. 
When the collapse center reaches a density that yields a 
Jeans length equal to twice the finest grid size, we terminate 
the simulation; typically the density dynamical range is 10 ^®. 
Periodic condition is employed throughout this study. 

We treat the two-phase gas as a single-component gas 
with an unusual equation of state. The index of adiabaticity 
7 approaches 0.7 in the interior of the molecular core and 1.1 
for the almost isothermal warm HII gas. For any gas dynam¬ 
ics simulation code written in a conservation form, the inter¬ 
nal energy density e and the pressure P are both needed as 
the dynamical variables. We model the equation of state as 
P — a{s)g{p), and the squared sound speed = 'y{p){P/p), 
where p the mass density, a{s) the entropy function, s is the 
entropy and the adiabatic index 7 (p) = ding{p)/d\np. To 
find the internal energy density, e, from the pressure P, the 
first law of thermodynamics demands. 


de 

P 


e + P 


dp — 0 , 


( 1 ) 


and the internal energy density e is thus related to the pres¬ 
sure P through 


e = P 


\j;-^ r#* 

Lff(p) Jo ® 


( 2 ) 


To model the change in the adiabatic index in the two- 
phase gas, we first let g{p) = {p/poY^'’\ and given the lim¬ 
iting values of 7 (p) in the core and in the warm gas, we 
model 


7(p) = 7 q — A7tanh 



(3) 


where the core density pc = po exp(1.0277) and the warm gas 
density pg = po exp(—1.07p). Here, po is the density at the 
transition boundary of the two-phase gas, and p represents 
the range of transition density in unit of po, chosen to be 
p = 3.3 to yield pdpg ~ 1000. Finally, we set 70 = 0.9 and 
A 7 = 0.2. (Note that 7 (p) is different from 7(p) due to an 
extra contribution from dy/dp.) 
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Figure 1. The density dependence of ln(e/a(s)) (red dashed line), 
ln(P/a(s)) (blue dashed line), 7 (p) (red solid line) and 7 (p) (blue 
solid line). 

Fig. 0 depicts various thermodynamics quantities, 
P/a{s),e/a{s),'y and 7 , as a function of density. Note that 
7 and 7 closely track each other, and the ratio between e 
and P increases with the density due to an increasing heat 
capacity. 

We start the simulation with the aforementioned 
uniform-density, quasi-spherical H 2 core of mass 1.5 Mq 
and temperature 20 K, which is in pressure balance with 
the surrounding warm gas of 10"^ K. (A ’’quasi-sphere” here 
is a sphere superposed with about 20% of higher multi-pole 
moments.) So long as the core size is smaller than the Jeans 
length, a uniform core is stable, and the core oscillates with 
large amplitudes since the gravitational force is initially un¬ 
balanced. The oscillation period is about 2 Myr. 

The core can quickly adjust to establish a pressure gra¬ 
dient within the initial 10 oscillation periods to balance the 
gravity. The central temperature of the core by now drops 
from 20 K to 12 K and the central density rises from 10^ 
cm“® to 5 X 10"^ cm“®. After that, the enhancement of the 
central density slows down while the core continues to os¬ 
cillate with finite amplitudes. When the central density be¬ 
comes sufficiently large after about 100 oscillation periods, 
the core becomes unstable and spontaneously collapses. The 
evolution of the central density in Fig. © shows the last tens 
of oscillations prior to the spontaneous collapse. 

The peculiar feature of a strong density gradient for 
driving the core instability will be discussed in the next 
section. For now we shall address the mechanism by which 
changes of the density gradient are possible. 

An isentropic equilibrium with spherical symmetry can 
only support potential flows, and it cannot produce vortex 
flows to yield material mixing. Hence the density gradient 
cannot change securely. However, an equilibrium sphere with 
an existing entropy gradient may generate vortical motion 
and can lead to material mixing. Rayleigh-Taylor instabil¬ 
ity is a familiar example. In our case, the system has an 
entropy gradient at the core boundary. But it is Rayleigh- 
Taylor stable since the cold matter is at the bottom of the 
gravitational potential instead of on the top. Nevertheless, 
we have a core on one hand possessing high-order multipole 
moments and on the other hand oscillating with large am¬ 
plitudes. Misalignment of the entropy and density gradients 



tvo 


Figure 2. Time evolution of the central densities of quasi¬ 
spheroid (thick dashed line) and ellipsoid (thin dotted line) nor¬ 
malized to 10‘^Pg, the density gradient A{t)/A(fi) evaluated by 
Eq. (|8ll (thick solid line), the first two terms (thin dashed line) 
and all terms (thick dash-dotted line) on the right of Eq. 0. 
The core oscillates and relaxs to a final configuration accessible 
to collapse. 

can develop as a result of nonlinear coupling of oscillations 
of different multi-pole moments. Therefore weak vortices are 
generated, thereby mixing the core over a long time. 

To demonstrate, we show a slice of the core through 
the core center in Fig. dla), and depict the two-dimensional 
vortical velocity normalized to the local sound speed at the 
moment of 72 oscillation periods. The vortical velocity field 
is computed by a projection operator I — fcfc acting on the 
Fourier component of the velocity field, where I and k are the 
identity tensor and the directional unit vector of wavenum¬ 
ber k, respectively. Sub-sonic weak vortex flows are clearly 
seen, where the Mach number is about few percent. On this 
slice, the flow exhibits a quadrupole pattern, for which a ra¬ 
dial inflow on one axis is accompanied by a radial outflow 
on another axis, and the flow has a negligible component in 
the direction of the third axis. Along the third axis, the flow 
velocity is also small. Thus, the flow is on a thick torus in 
three-dimensional space, and it has a dominant quadrapole 
velocity pattern. The direction of the vortex flow can actu¬ 
ally reverse and it oscillates quasi-periodically with frequen¬ 
cies comparable to the density oscillation. 

To investigate density and entropy, the respective spher¬ 
ical shell averages are subtracted from these two quantities. 
The contours of the instantaneous high-order multipole mo¬ 
ments of density and entropy are plotted in Fig. (|3lb) on the 
same slice as Fig. dSla). The density and entropy gradients 
are clear in regions where vortex flows are strong, and the 
misalignment between these two gradients are the source of 
oscillating vortex flows. 

In Fig.(|31c), we plot the spherical shell average profiles 
of density and entropy. The scatters around the average pro¬ 
files indicate the magnitudes of high-order multipole mo¬ 
ments. The initial sharp density and entropy boundary is 
now smeared over an extended regions after 72 oscillation 
periods due to the material mixing, which yields a more 
amenable configuration for a molecular core. Therefore we 
shall take the relaxed configuration of Fig.([3lc) as the initial 
condition for the study to follow. 
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dient A, the first two terms and all terms combined on the 
right of Eq. 0 in Fig. ©• As the core contracts, the aver¬ 
age density gradient A increases and the other two curves 
also increase. All three curves oscillate in phase, indicative 
of that the restoring force associated with the increase of 
/^-average pressure is capable of resisting the enhancement 
of self-gravity and produces stable nonlinear oscillations. 

The above analysis of the hrst-order virial equation for 
nonlinear stability is at best ad hoc, largely relying on the 
simulation result. As the stability of such a soft core is new 
and counter-intuitive, we now analyze this problem in a more 
rigorous manner. Despite the core configuration in the sim¬ 
ulation is away from the dynamical equilibrium, the under¬ 
lying linear stability of the core must be obeyed to account 
for the simulation result. We thus assume the existence of a 
quasi-static equilibrium and perform a linear stability analy¬ 
sis about this equilibrium. To simplify the analysis, we model 
the dynamical equilibrium to consist of a cold spherical core 
possessing a density gradient and a uniform warm diffuse 
gas that is separated from the core by a sharp boundary. 

In a dynamical equilibrium, the first-order virial is iden¬ 
tically zero. Variational principle deals with the second-order 
virial, and the second-order potential energy S^W of linear 
perturbations can be written as, 

= j dV[7Po(V ■ ■ Vpo)(4 ■ VM 

+ 2(^ ■ VPo)V • ^ -b (po^ • V<5<^)] (5) 


Figure 3. (a): A slice of the vortical velocity field normalized 
to the local sound speed. This slice cuts across the center of the 
core and has a size of 6Rc. The maximum Mach number is about 
4%. (b): The contour of the high-order multipole moments of 
density (dark curve) and entropy function P/g{p) (light curve) 
on the same slice of (a). The density and the entropy function are 
normalized to the background values in the HII gas. The density 
and entropy gradients are the source of oscillating vortex fiow and 
are co-spatial with the vortex flow, (c): The shell-averaged radial 
profiles of density (solid line) and the entropy function P/g(p) 
(dashed line). The normalizations are the same as (b). The error 
bars indicate the magnitudes of high-order multipole moments. 


3 STABILITY OF A SOFT CORE 

Despite the very soft equation of state, the core is found to 
be nonlinearly stable. This counter-intuitive core stability 
can be seen from the virial equation, 

=3(P)-3P(P) + (D). (4) 

where the volume average is defined as (/) = 

fd^r/{4:TiR^/3) integrated up to a fixed radius R, which 
we call P-average, ^r{r,t) is the radial displacement, Q the 
potential energy density, and P(P) is the pressure evaluated 
at radius R. For convenience, we choose R to be located 
immediately outside the core radius. With 7 < 1, the pres¬ 
sure within the core must decrease upon contraction. How¬ 
ever, another contribution to the P-average pressure from 
the warm gas does the opposite. It increases the P-average 
pressure when the core contracts since the volume occupied 
by the warm gases within P increases from zero to a pos¬ 
itive value. To demonstrate the stability, we compute from 
the simulation the time evolution of the average density gra- 



dP-^[(^-VPo) + 7Po(V-0], 


where the surface integral is evaluated at P, just outside 
the core radius Pc. The lowest-energy mode must avoid ex¬ 
citing the positive-energy sound wave in the high tempera¬ 
ture region, and hence the displacement is incompressible, 
V ■ ^ ^ 0. Moreover, as the density is low and the gravity 
negligible, the warm gas has a vanishingly small equilibrium 
pressure and density gradients, i.e., VPq —>■ 0 and Vpo —^ 0 . 
The incompressibility of the warm gas also yields a van¬ 
ishing small density perturbation, thus a vanishingly small 
gravitational perturbation, S<j) 0. Finally, the surface in¬ 
tegral is also small due again to incompressibility and to the 
small pressure gradient just outside the core. The second- 
order potential energy 5^W is thus contributed entirely by 
the region r ^ P. 

Like a pulsating star, the most unstable mode involves 
only the radial displacement, ^ = ^rr, so that it can tap the 
gravitational energy. Eq. 0 then becomes 


rW =471 




d{Y^r) 

dr 


2 c 2 dpQ d(j>o 


+ 


dPo d(r^Cr) 
dr dr 


- r^5c47rGpoj. 


( 6 ) 


Note that dpo/dr in the second term contains a Dirac-5 func¬ 
tion at the core boundary. This term is positive and results 
in the stabilizing effect mentioned earlier, i.e., an increase of 
average pressure within a fixed P contributed by the warm 
gas to stabilize the perturbations. 

We next parametrize the equilibrium density as poix) = 
po(l)(a; -I- A{1 — x)), where x = r/R and A > 0, i.e., a linear 
model from which the force balanced pressure can be calcu¬ 
lated. It yields Po{x) = po(l)2o(l) -|- 7 rGpo(l)^P^[(l/ 4 ) -|- 
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(5/18)yl + (5/36)A2 - {2/2,)A^x^ + (7/9)A(A - l)x® - 
(l/4)(yl — Here, A can be regarded as the average 

density gradient of the core. 

Finally we adopt a trial function — (3/4)a;'^ — 

I3{x^ — (3/5)a;®). Here /3 is an optimization parameter, where 
5^W is to be minimized with respect to /?. An optimal choice 
of the trial function can provide an accurate estimate for the 
marginal stability condition. This trial displacement gives 
uniform compression at a: = 0 and becomes incompressible 
at a: = 1. Demanding the radial displacement and its ra¬ 
dial derivative to be continuous across the boundary, we let 
P^r = 1/4—(2/3/5) outside the core, which is incompressible 
as desired. 

Substituting the above into Eq. © and setting the first 
/3-derivative of 5^W equal to zero, we find that the optimal 
/3 = Pop- For all practical purposes when Pop is on the order 
of unity, S^W is insensitive to p. Hence for the time being 
we set /3 = 0 to simplify the algebra and the presentation. 
Full expressions due to a finite /3 will be given in Appendix. 

Eq. ® becomes 


Pw = lMHi 6770(1) , 

Rl I 5 


nGR^Po(l) 


420 


/(111 + 100.1 + 29^) - + + j 

17) 


where we have let R = Ro- 

On the other hand, the core mass, denoted as Mo, is a 
conserved quantity and related to the core radius R^ by the 
following relation. 


_ OlVlc - , . 

" [^po(l)(3 + A)J • ^ ^ 

We define "Fo)!) = —GMo/Rc- Substituting the Eq. I® into 
Eq. (jT]), the unstable condition, where S^W is negative, be¬ 
comes. 


Toil) 


[(f-29,).7. 

- IOO 7 ) A -b (120 - III 7 )]. 


■ 168(3-bA)7 
'695 


(9) 


Note that core relaxation increases the average density 
gradient A. As the right-hand side is an increasing function 
of positive A for 7 less approximately than unity, the core, at 
a fixed ro(l)/|'Fo(l)|, therefore becomes more unstable when 
the density gradient A increases. The ratio To(l)/|4>o(l)| is 
measured to be 1.28 from the simulation at the moment of 
spontaneous core collapse, and the marginally stable value 
of A is thus predicted to be 8.9 for 7 = 0.7 according to Eq. 
® . (For a finite /3 the marginally stable A is predicted to be 
8.75.) This is to be compared with the measured A = 7 from 
the simulation at core collapse. Given the crude linear model 
with a sharp boundary for the linear stability analysis, the 
prediction is in fair consistency with the simulation result. 

The evolution of A also shown in Fig. m reveals that 
the core configuration begins with A = 5.6 and evolves to¬ 
ward collapse when A = 7 in 28 initial oscillation periods. 
Therefore, the small density gradient explains the stability 
of soft core and may explain the stability of Barnard 68 . 

The measured A is found by the following procedures. 
The core mass Mo is first evaluated by the configuration 
before core relaxation where the core radius Ro is chosen to 


be the radius of a uniform core. Since Mo is conserved, the 
core radius Rc can be determined by calculating the enclose 
mass at any moment. The density at the boundary po(l) can 
be measured once Rc is known. Finally, the density gradient 
A is determined by using Eq. ®. 


4 CORE COLLAPSE TRIGGERED BY 

SMALL-AMPLITUDE BAGKGROUND 

SOUND WAVES 

Two different cases are to be investigated. First, a stand¬ 
ing wave acts on a quasi-spherical core studied in the last 
section; this case provides a baseline stability criterion. Sec¬ 
ond, a Gaussian pulse impinges on an elliptical core; this 
case addresses a more general situation. 

Since the core is oscillating, it is natural to examine the 
interaction efficiency as a function of wave frequency. Here 
the core is located at the pressure node of the standing wave. 
When wave frequency coincides with the core frequency, the 
interaction is expected to be most efficient, thus requiring 
the lowest wave amplitude to drive the collapse. We define 
the threshold wave amplitude as the one that can marginally 
compress the core to collapse. The marginally driven case 
is one that the core peak density rises monotonically, but 
hesitantly, toward collapse in one to two oscillation period. 

We find that the resonant interaction is indeed the most 
efhcient; however the quality factor of resonance is rather 
low, as shown in Fig. ®. Here, we plot the threshold wave 
amplitude as a function of the wave frequency and the rel¬ 
ative phase between core and wave oscillations. Since col¬ 
lapse occurs in one core oscillation period, one therefore 
expects the relative phase to be also a crucial parameter. 
Fig. (0 contains 6 groups of data points. Data in every 
group have the same wave frequency, and six frequencies 
V = 0.5, 1,1.5, 2, 2.5 and 3 times of the initial core oscilla¬ 
tion frequency vq of Fig. m are investigated. Data in the 
same group represent different relative phases, and data of 
the same symbol in different groups represent the same rel¬ 
ative phase but different wave frequencies. 

Since wave and core frequencies are different, the rela¬ 
tive phase needs to be defined carefully. Here, it is the wave 
phase when the core oscillation phase is most rarified, or 
vice versa. The zero relative phase p = 0 corresponds to the 
situation that the most rarified phases of the two oscilla¬ 
tions coincide. A non-zero p refers to the phase of the wave 
oscillation or the core oscillation whichever is faster. 

It is seen that (/ = 0 is almost always the most efficient 
phase for all frequencies to drive core collapse. When (/ = 0 
and u = l.Or'o, it requires the lowest wave amplitude to drive 
core collapse. When </> = 0 but z/ yf z/q, the mismatch in os¬ 
cillation timing requires a more vigorous sound amplitude 
to drive core collapse. However, it turns out that not much 
higher amplitude is needed to drive the off-resonance col¬ 
lapse, a clear indication of broadband resonance. The 3dB 
bandwidth of the threshold wave power is about a factor 
6 , from zz/z/Q = 0.5 to 3, with an allowable range of phase 
mismatch about ±60 degrees. 

We next consider an ellipsoidal core impinged by a 
Gaussian sound pulse to test whether the above result holds 
in a more general case. The ellipsoid is parameterized by the 
lengths of the three principle axes, which are respectively 
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Figure 4. Frequency/phase response of the threshold wave am¬ 
plitude {6pg/pg) for triggering collapse of a quasi-spheroid core 
by a sinusoidal standing wave. The six groups of symbols cor¬ 
respond to ly/i^o = 0.5,1,1.5, 2, 2.5, 3, where ly and i^o the 
wave and core frequencies, respectively. Different symbols in each 
group represent different initial relative phases between the two 
oscillations. 


l.SRc, Rc and Rc/l.S, and obtained by uniform volume- 
preseving distortion of the above quasi-spherical core, so 
that the core mass is preserved. The ellipsoidal core is also 
stable and undergoes similar relaxation oscillations as the 
quasi-spherical core, albeit the stable oscillation ends earlier 
and the amplitude is twice greater. Fig. (H) also shows the 
peak density evolution of the ellipsoidal core as it relaxes to 
collapse. 

A Gaussian sound pulse of half width cr contains a wide 
range of wavelengths, with a characteristic wavelength about 
twice the full-width, 4(j, and the characteristic frequency 
v = Cs,gl^<y, where Cs^g is the sound speed of the warm 
gas. The relative phase can be defined to be proportional 
to the ratio of the distance of the pulse waist from the core 
center to the sound traveling distance I/o(= C's,g/r'o) in one 
core oscillation period. Hence, we use these two parameters, 
I/o/4(t and L/Lq, to replace the wave frequency and the 
relative phase of the previous case to search for the threshold 
pulse amplitudes. To be exact, we define the phase as cf) = 
27r(0.34— (L/Lq)) so that the interaction is most effective at 
(j) = 0, and the pulse waist is at 1.276 (j from the pulse peak 
toward the core. The two numerical values, 0.34 and 1.276, 
are determined after the simulation results are analyzed. In 
this test, the Gaussian planar sound pulse propagates along 
the major axis of the ellipsoid. 

Fig. ([5)) depicts the result. It follows the general trend 
as that obtained in the previous case, except for about three 
times higher threshold wave amplitudes. The ellipsoid cre¬ 
ates phase incoherence for the incoming pulse; even when 
the pulse can reach the core center, the distorted wave front 
arrives asynchronously, thus requiring a higher wave ampli¬ 
tude to trigger the collapse. Nevertheless, the dependence on 
the relative phase between core and wave oscillations is now 
weaker than the previous case. When the wave frequency 
meets the broadband resonance condition, the triggered col¬ 
lapse can take place with more than 50% efficiency. 

We also show in Fig. m the radial profiles of density, 
mass flux and Mach number of the infall for the ellipsoidal 
core at three instants immediately before the collapse. The 



v/vq 

Figure 5. Same as Fig. (l4]l except for the target being an ellip¬ 
soidal core and the driver being a Gaussian pulse. The charac¬ 
teristic frequency of the Gaussian pulse and the arrival phase are 
defined in the text. 



r/Rc 


Figure 6. A time sequence of shell-averaged radial profiles of den¬ 
sity (diamond), infall Mach number (square) and radial mass flux 
(dot) near the instant of collapse. The time sequence is marked 
from light symbols to dark symbols. The mass flux is normalized 
to WOCg ,,/G for the clarity of display. The error bars represent 
deviations from spherical symmetry. 


density profile fits approximately by a power law with a 
— 1.80 logarithmic slope, and the mass flux profile is almost 
flat with a logarithmic slope —0.10. The Mach number ex¬ 
ceeds unity at the radius when the mass flux begins to fol¬ 
low the power law, indicating that the power-law mass flux 
is established in the pressure-free regime. By contrast, the 
density follows the power-law behavior at a much greater 
radius when the core pressure is still important. The mass 
flux, or m, is about C Y/G, com parable to the prediction 
of the inside-out model IlShull 19771 ). where Cs.c is the sound 
speed at the core boundary. However, the ev olution of the 
core is consistent with the o utside-in collapse (jLarso ni l 19691 : 
IPenstorJl 19691 : iHunterll 1977h . Although profiles for the spon¬ 
taneously collapsed core described earlier are not shown here 
but we find that their general behavior is similar to the 
driven case, but the mass flux is a factor two greater. 
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5 MULTI-CORE COLLAPSE TRIGGERED BY 

A SINGLE WAVE 

It is interesting to consider the possibility that the passage of 
a single sound pulse may trigger collapses of multiple cores; 
this may yield propagating star formation. Since the favor¬ 
able sound wave has a wavelength few tens greater than the 
core and since the wave front can extend transversely over a 
large distance, the wave can in principle interact with a num¬ 
ber of cores at different locations and different times. We let 
a single Gaussian pulse interact with five nearby cores which 
are placed in a region of size about Icr. The five cores are 
identical ellipsoids considered in Sec. The core oscilla¬ 
tion frequencies are near the optimal for collapse but the 
phases are random. 

The simulation result is shown in Fig. 0 . As a pulse 
slightly exceeding the marginal strength arrives, three cores 
that are out of the favorable phase range do not collapse, but 
two cores that do collapse. The two cores collapse in different 
manners. One (a) that is near the optimal phase collapses 
promptly, and the other (b) that somewhat deviates from 
the optimal phase hesitates for a while before it decides to 
collapse. Whether these cores collapse or not agrees with the 
criterion given in Sec. @ for an isolated core as if no other 
cores were to exist. 

The result may be somewhat puzzling, since back- 
reactions from the cluster of impinged cores seem not to 
matter. We find it can be understood as follows. The sound 
inside the core travels much slower than the Gaussian pulse. 
Hence any back-reaction wave emanated from the interior 
of the impinged core will take a finite time to get out of the 
core, and by the time the perturbation reaches the neigh¬ 
boring cores, the Gaussian wave has already acted on these 
cores. Hence, the individual core collapse agrees with the 
criterion given in Sec. 


6 CONCLUSIONS 

We consider the linear and nonlinear stabilities of a two- 
phase gas in the low-mass star forming regions, where a 
dense core of adiabatic index 7 = 0.7 is surrounded by an al¬ 
most isothermal warm diffuse gas of 7 = 1.1 and in between 
is a blanket of transition layer. 

We show that a soft dense core can be stably oscillating 
as a result of a small average density gradient. A molecular 
core of 10 K temperature, 0.2 pc size and 1.5 Mq mass yields 
a typical oscillation period about 2 Myr. However the core 
will relax and slowly evolves to spontaneous collapse in tens 
of Myr if left unperturbed. This relaxation of dense core is 
found possible due to the presence of a vortex flow, which is 
generated by the very existence of entropy gradient around 
the soft dense core. 

We have also demonstrated a triggered collapse mech¬ 
anism of a stably oscillating soft core by means of a small- 
amplitude background wave. This mechanism involves reso¬ 
nant interactions, but the resonance is broad-band, indica¬ 
tive of low-quality factor resonance. The conditions for the 
triggered collapse are that wave frequencies are comparable 
to and up to a factor a few higher than the core oscillation 
frequency, and that the encountering phase should be within 
about ±90 degrees of the optimal phase. The latter condi- 


(a) 

(b) 

(a) 

(b) 

(a) 

(b) 

• 

• 0 

. 

(b) 


Figure 7. Four sequential snapshots of a Gaussian sound wave 
propagating from left to right interacting multiple ellipsoidal 
cores. It is noted that two ((a) and (b)), out of five, cores collapse 
during the wave passage. The five cores have identical oscillation 
frequencies but differ by the wave encounter phases. 


tion means that this mechanism has about 50% efficiency 
even when the resonance condition is met. 

The triggered collapse mechanism however requires long 
waves. Sound waves in the warm HII region of 10, 000 K 
must have a wavelength ranging from 6 to 36 pc to produce 
efficient interactions. These small-amplitude long waves may 
be originated from active sources located tens of pc, or even 
one hundred pc, away from the star forming regions through 
the push of ionization fronts and outflows that varies on the 
time scale of 0.7 — 4.5 x 10® years. 

The passage of the sound wave can trigger several col¬ 
lapsing cores at different locations and different times. Thus, 
stars may form with tight spatial-temporal correlation trac¬ 
ing the propagation of the wave front, yielding propagating 
star formation for cores within a range of masses. If the tem¬ 
perature of different cores remains roughly about 10 K, the 
oscillation period will scale linearly with the core size, and if 
these cores are not much below the Jeans instability bound¬ 
ary, they will have masses scaled linearly with the core size, 
thus the oscillation period. A single wave passage can there¬ 
fore trigger collapse of cores within a mass range slightly 
short of one order of magnitude. 

It follows that young stellar objects (YSO) may have 
spatial-temporal correlation. This prediction can be tested 
by measuring the locations, ages and spectral types of YSO 
in a region of 100 pc. Over this region, the sound transit time 
is on the order of 10^ years. This time scale coincides the life 
time of YSO before the YSO enters the main sequence, dur¬ 
ing which the time stamps of various evolutionary stages are 
distinct, permitting the propagating star formation scenario 
to be testable. 

In this work, we do not consider the effects of rota- 
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tion and magnetic fields. Our preliminary results show that 
weak rotation of the core does not affect the triggered col¬ 
lapse mechanism reported in this study; it only changes the 
central configuration of the collapsed core into two collapsed 
centers. Magnetic fields may further complicate the central 
configuration. We will report these effects in a series of fu¬ 
ture works. 
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APPENDIX 

In Sec. (|3]), we have introduced an optimization parame¬ 
ter (3 to minimize 5^W, which yields the optimized f3op = 
G{A, 7 , a)/2F{A, 7 , a), where G and F would be defined in 
the next paragraph. We find that /? = 0 offers a good approx¬ 
imation to the solution to this problem, and therefore take 
/? = 0 to simplify the presentation in the main text. Here we 
give a full expression of S^W for the finite optimized value 
of Pop- 


The full expression of S^W, including (3, can be ex¬ 
pressed as 

TV^Gplil) 




Rc 


■[F{A,j,a){P-popf+ 

AF{A,y,a)H{A,y,a)-G\A,'y,a) 
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where a = ro(l)/|$o(l)| of Eq. (|9l), 
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(13) 


For S^W to have a minimum with respect to f3, we need 
F > 0, and for the marginal stability, we let AFFI — G^ = 0. 
The condition AFFl — G^ — 0 is found to yield E > 0 for all 
positive A, 7 and a, hence justifying the /3 parametrization. 

Fig. ([ 8 ]) plots the family of AFH — = 0 as a function 

of 7 and a with A as the family parameter. Above each 
curve, 4FFF — G^ > 0, representing the stable region. As one 
can see, the family of curves asymtote to a straight line for 
which A —>■ 00 and 7 has a limiting value 1.309. 

This stability criterion indicates that when 7 increases, 
harder cores become more stable and one needs a larger 
A for the core to be unstable. However the stability cri¬ 
terion can become invalid when 7 = 1.309, where it pre¬ 
dicts that the marginal stable A is infinity, and when 7 is 
greater than this value, the core becomes completely stable. 
Clearly the result cannot be correct. The number 7 = 1.309 
should not be taken literally, but the tendency of a hard 
core for the above stability analysis to be invalid may be 
understood as follows. For a hard core, the equilibrium den¬ 
sity p{r) tends to cross zero at a finite radius, the profile 
is convex and hence it should be better approximated by a 
quadratic model, e.g., po{x) ~ 1 -I- (A — 1)(1 — x^). On the 
other hand for a soft core, the equilibrium density tends to 
extend to infinity, the profile is concave in most part of the 
mass and may be captured by the linear model. For exam¬ 
ple only when 7 ^ 1.2, the Lane-Emd en equation permits 
a pow er-law density p(r) ~ i'll dBinnev fc Tremainel 

l2008ll . the profile is concave and we think that Eq.® can 
be approximately valid. Extension of the above approximate 
treatment of linear stability analysis to a hard core is a in¬ 
teresting problem in other applications and worths further 
investigations. 














Sound-Triggered Collapse of Stably Oscillating Low-Mass Cores in a Two-Phase Interstellar Medium 9 



Figure 8 . Curves of AFH — = 0 as a function of 7 and a. 

Different curves represent different values of A. The stable region 
lies above each curve. For A —>• 00 , 7 has an upper bound value 
1.309. 















